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Abstract 

Large-eddy simulations of incompressible Newtonian fluid flows with approximate 
deconvolution models based on the van Cittert method are reported. The Legen- 
dre spectral element method is used for the spatial discretization to solve the fil- 
tered Navier-Stokes equations. A novel variant of approximate deconvolution mod- 
els blended with a mixed scale model using a dynamic evaluation of the subgrid- 
viscosity constant is proposed. This model is validated by comparing the large-eddy 
simulation with the direct numerical simulation of the flow in a lid-driven cubical 
cavity, performed at a Reynolds number of 12'000. Subgrid modeling in the case 
of a flow with coexisting laminar, transitional and turbulent zones such as the lid- 
driven cubical cavity flow represents a challenging problem. Moreover, the coupling 
with the spectral element method having very low numerical dissipation and dis- 
persion builds a well suited framework to analyze the efficiency of a subgrid model. 
First- and second-order statistics obtained using this new model are showing very 
good agreement with the direct numerical simulation. Filtering operations rely on 
an invertible filter applied in a modal basis and preserving the C*^-continuity across 
elements. No clipping on dynamic parameters was needed to preserve numerical 
stability. 
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1 Introduction 



Large-eddy simulation (LES) represents a way of reducing the number of de- 
grees of freedom of the simulation with respect to the requirements of the direct 
numerical simulation (DNS). This is done by calculating only low- frequency 
modes in space and modeling high-frequency ones, the scale separation be- 
ing performed by filtering in space the Navier-Stokes equations. Large-scale 
structures are obtained by the computed flow dynamics while the behavior 
of subgrid scales and their interaction with large eddies are modeled by the 
additional term in the LES governing equations resulting from filtering the 
Navier-Stokes equations. The expression of the additional term as a function 
of the resolved field is referred to as subgrid modeling. 

Approximate deconvolution models (ADM) constitute a particular family of 
subgrid models. They rely on the attempt to recover, at least partially, the 
original unfiltered fields, up to the grid level, by inverting the filtering operator 
applied to the Navier-Stokes equations. The focus here is on the approximate 
iterative method introduced by Stolz and Adams [1] which is based on the 
van Cittert procedure. This method was subsequently applied to incompress- 
ible wall-bounded fiows [2], to compressible fiows and to shock-boundary layer 
interaction [3] using a new variant ADM-RT, blending ADM with a relax- 
ation term (RT) increasing the dissipative character of the model. Transitional 
flows were also investigated by Schlatter et al. [4]. Over the past flve years, 
ADM spread over various fields of application. Gullbrand and Chow studied 
the effect of explicit filtering in the case of channel flow [5]. ADM were also 
more recently applied to the LES of a rectangular jet and to computational 
aero- acoustics by Rembold and Kleiser [6] . Particle-laden turbulent flows were 
investigated in the ADM framework by Shotorban and Mashayek [7]. Prom 
the numerical viewpoint, Schlatter et al. [4] used a parallel implementation 
of a mixed Fourier-Chebyshev spectral method. These models were also im- 
plemented in a flnite volume framework in the semi-industrial code NSMB, 
Navier-Stokes Multi-Block, by von Kaenel et al. who apphed it to shock- 
boundary layer interaction and channel flow in [8,9]. To our knowledge, the 
only implementation based on the spectral element method (SEM) is due to 
Iliescu and Fischer [10] who used ADM based on the rational LES model 
(RLES) instead of the van Cittert one. More recently, Pruett et al. proposed 
a temporal ADM for LES [11] and a stability analysis of the LES- ADM equa- 
tions was performed by Dunca and Epshteyn [12]. 
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LES of Newtonian incompressible fluid flows with ADM based on the van Cit- 
tcrt method using Legendre-SEM as spatial discretization to solve the filtered 
Navier-Stokes equations are envisaged for the first time in this paper. Follow- 
ing the idea of Winckelmans et al. [13] who coupled the ADM based on the van 
Cittert method and the Smagorinsky model [14] , and GuUbrand and Chow [5] 
who proposed a dynamic version of the previous model, a new variant which 
blends ADM and the mixed scale model introduced by Sagaut [15] with a dy- 
namic evaluation of the subgrid- viscosity constant based on a Germano-Lilly 
type of procedure [16, 17] is proposed. 

A specific filtering operation adapted to SEM and preserving continuity across 
elements is applied in a modal basis which was proposed in the p-version of 
finite elements and first used by Boyd [18] as a filtering technique. Depending 
on the transfer function, this filter is not projective and is therefore invertible, 
this property being essential for the deconvolution procedure. 

A DNS of the fiow in a lid-driven cubical cavity performed at Reynolds num- 
ber of 12'000 with a Chebyshev collocation method due to Leriche and Gavri- 
lakis [19] is taken as the reference solution to validate the new model. Sub- 
grid modeling in the case of a fiow with coexisting laminar, transitional and 
turbulent zones such as the lid-driven cubical cavity fiow represents a chal- 
lenging problem. As the fiow is confined and recirculating, any under- or over- 
dissipative character of the subgrid model can be clearly identified. Moreover, 
the very low dissipation and dispersion induced by SEM allows a pertinent 
analysis of the energetic action induced by any subgrid model, which is not 
feasible in the framework of low-order numerical methods. The coupling of 
the lid-driven cubical cavity fiow problem with the SEM builds therefore a 
well suited framework to analyze the accuracy of the newly defined subgrid 
model. Bouffanais et al. in [20, 21] have performed LES of the fiow in a lid- 
driven cubical cavity at a Reynolds number of 12'000 using the same physical 
parameters as the DNS from Leriche and Gavrilakis [19]. The numerical frame- 
work of [20,21] is the same as the one used in the present article. Standard 
subgrid models were used in [20,21]: dynamic Smagorinsky [16,17] or dynamic 
mixed models [22]. 

The paper is organized as follows. In Section 2, the filtered Navier-Stokes equa- 
tions are given, followed by a brief description of the space-time discretization 
using the spectral element method. The subgrid modeling is dealt in details 
in Section 3 and the numerical filters are described in Section 4. LES of the 
fiow in the lid-driven cubical cavity, based on the subgrid models introduced 
in Section 3, is presented and thoroughly analyzed in Section 5. Finally, in 
Section 6 we present the conclusions. 
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2 Governing equations and numerical method 



2.1 Governing equations 



In the case of isothermal flows of Newtonian incompressible fluids, the LES 
governing equations for the flltered quantities denoted by an overbar, obtained 
by applying a convolution fllter Q-k to the Navier-Stokes equations, read 

dui d dp d f dui duj\ dry 

dt dxj * ^ dxi dxj \ dxj dxi J dxj ' 

the flltered velocity fleld u = ^ * u satisfying the divergence-free condition ([2]) 
through the flltered reduced pressure fleld p. The components of the subgrid 
tensor r are given by 

"Tij ~ UiUj UiUj, (3) 

and u is the kinematic viscosity. The closure of the flltered momentum equation 
([T]) requires r to be expressed in terms of the flltered fleld which reflects the 
subgrid scales modeling and the interaction among all space scales of the 
solution. 



2.2 Space discretization 



The numerical method treats Eqs. ([I])-([2]) within the weak Galerkin formu- 
lation framework. The SEM consists in dividing the computational domain 
into a given number of spectral elements. In each spectral element, the veloc- 
ity and pressure flelds are approximated using Lagrange-Legendre polynomial 
interpolants. The reader is referred to the monograph by Deville et al. [23] for 
full details. The velocity and pressure are expressed in the Pp — Pp_2 func- 
tional spaces where Pp is the set of polynomials of degree lower than p in each 
space direction. This spectral element method avoids the presence of spurious 
pressure modes as it was proved by Maday and Patera [24,25]. The quadra- 
ture rules are based on a Gauss-Lobatto-Legendre (GLL) grid for the velocity 
nodes and a Gauss-Legendre grid (GL) for the pressure nodes. 

Borrowing the notation from Deville et al. [23], the semi-discrete flltered 
Navier-Stokes equations resulting from space discretization are 

M^ + Cu + z/Au-D^P + Dt = 0, (4) 

-Du = 0. (5) 
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The diagonal mass matrix M is composed of three blocks, namely the mass 
matrices M. The global vector u contains all the nodal velocity components 
while p is made of all nodal pressures. The matrices A, D^, D are the discrete 
Laplacian, gradient and divergence operators, respectively. The matrix oper- 
ator C represents the action of the nonlinear term written in convective form 
u- V, on the velocity field and depends on u itself. The semi-discrete equations 
constitute a set of nonlinear ordinary differential equations (jlj) subject to the 
weak incompressibility condition ([5]). 

2.3 Time integration 

Standard time integrators in the SEM framework handle the viscous linear 
term and the pressure implicitly by a backward differentiation formula of 
order 2 (BDF2) to avoid stability restrictions such that vAt < Cjjy"., while 
all nonlinearities, including the discretized subgrid term — Dr, are computed 
explicitly, e.g. by a second order extrapolation method (EX2), under the CFL 
restriction 

MmaxAt < C/pl (6) 

The implicit part is solved by a generalized block LU decomposition with 
a pressure correction algorithm [23,26,27]. The overall order-in-time of the 
afore-presented numerical method is two. 



3 Subgrid modeling 

3.1 General considerations 

The problem of subgrid modeling consists in taking into account the interac- 
tion between resolved and subgrid scales which is represented by the subgrid 
term V ■ r in the filtered momentum equation ([T]). 

Following the terminology introduced by Sagaut [28] , two modeling strategies 
are defined. A first group of models, called structural, aims at making the best 
approximation of the tensor r by reconstructing it formally from the resolved 
field u. The closure consists in finding a relation such that 

r- = C.(u), (7) 

where the upper index 'm' distinguishes the modeled from the exact subgrid 
tensor. This group of models does not require any foreknowledge about the 
nature of the interactions between resolved and subgrid scales. The second 
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group, called functional, consists in modeling the action of subgrid scales on 
the resolved field u using physical concepts and not at approximating the 
subgrid tensor r itself, even if a subgrid tensor is explicitly constructed as 
for subgrid-viscosity models. Most of these models assume that the action of 
subgrid scales on resolved ones is essentially energetic, so that the balance 
of energy transfers between both scales categories is sufficient to describe the 
interaction. 

The focus hereafter is on ADM which attempts to recover, at least partially, the 
original unfiltered fields, up to the grid level, by inverting the filtering operator 
applied to the Navier-Stokes equations. Following the idea of Winckelmans 
et al. [13] who coupled the ADM based on the van Cittert method and the 
Smagorinsky model [14] , and Gullbrand and Chow [5] who proposed a dynamic 
version of the previous model, a new variant blending ADM and the dynamic 
mixed scale model introduced by Sagaut [15] is proposed. 



3.2 Approximate deconvolution model 



The deconvolution approach aims at reconstructing the unfiltered fields from 
the filtered ones. The subgrid modes are not modeled but reconstructed us- 
ing an ad hoc mathematical procedure which falls in the structural modeling 
category. Writing formally the Navier-Stokes momentum equation ([T]) as 

^ + f (u) = 0, (8) 
the evolution equation of the filtered quantities becomes 

^ + f(u) = [f,^^](u), (9) 

where the convolution filter Qic = {Co'P)-k embodies the LES filter C-k and the 
projective grid filter V-k [5,13], the latter being therefore implicitly accounted 
for in the sequel. It is important to note that the LES filter and the grid filter 
do not commute since the effect of the SEM discretization is not a spectral 
cutoff filter, unlike the case of spectral methods as reported by Gullbrand and 
Chow [5]. The subgrid commutator reads then 



[f, Gk]{u) = i{Gku)-Gk f (u) = f (u) - f (u), (10) 
which is strictly equivalent to Eq. ([1]) given 

[f,^^](u) = -V-T. (11) 

The exact subgrid contribution appears as a function of the non-filtered field, 
which is not computed when performing a LES. This field being unknown, the 
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idea is to approximate it using the following deconvolution procedure 



u ~ u* = Qat 7^ u = (Qat o ^) 7^ u = {Qn o CoV)i^u = {Qn o C)-ku, (12) 

where u = V-ku is the grid- filtered velocity. The operator Qat* is an A^th-order 
approximation of the inverse of the filter C-k, since the grid filter is projective 
and therefore not invert ible, such that 

(Q^o/:)=J + 0(A^), (13) 

with Xt^t the identity filtering operator and A the filter cutoff length associated 
to Q-k. Stolz and Adams proposed in [1] an iterative deconvolution procedure 
based on the van Cittert method. If the filter C-k has an inverse, it can be 
computed using the truncated van Cittert expansion series 

C-'c^Qr, = j:{X-Cy, (14) 

which is known to be convergent if 

||J-£||<1. (15) 

The deconvolution error induced by the approximation f|T4l) can be represented 
by a filter Hn-^ defined by 

nN=I-QNoC. (16) 
The subgrid term is then approximated as 

[i,gMH ^ [i,g^]iQN^u) = [f,6;^](u*), (17) 

and the model resulting from this approach is obtained by introducing Eq. 
(ITTl) into the filtered Navier-Stokes momentum equation ([9]) 

^ + f(u) = [f,^^](u*). (18) 

Using once more approximation ( JT2l) in Eq. (fTSll implies f(u) = f(u*) and 
leads to the formulation commonly used with ADM 

^ + ^^f(u*)=0. (19) 

It is noteworthy that this latter formulation introduces the deconvolution er- 
ror and the error related to the non-inversion of Vk into the nonlinear ad- 
vection term, thereby breaking the Galilean invariance [29]. Furthermore, the 
expression of the subgrid tensor of Bardina's scale similarity model [30] is not 
recovered from the deconvoluted formulation (HM when Qn = 2, which is 
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again due to the difference between the filtered and the deconvoluted veloci- 
ties. Based on the previous comments, the filtered formulation (fT8|) appears to 
be the most general and therefore, all LES presented in the sequel rely on it. 
No numerical instabilities were observed using the formulation fllSp associated 
with our explicit treatment of the nonlinear terms, see Sec. 12.31 

3.3 Coupling with a dynamic mixed scale model 

Coupling ADM with a subgrid-viscosity model can be formally achieved by 
adding a source term s(u) to the right-hand side of Eq. (fT8l) 

^ + f(u) = [f,^^](u*) + s(n), (20) 

where s(u) is expressed in terms of the filtered rate-of- strain tensor S by 

s(u) = V ■ (z/sgs(Vu + Vu^)) = V ■ (2z/sgsS), (21) 

the superscript 'T' denoting the transpose operation and v^gs the subgrid vis- 
cosity. For such functional models, only the deviatoric part of the subgrid 
stress is modeled. On the other hand, the ADM part [f, ^★](u*) includes both 
isotropic and deviatoric parts. Using such subgrid-viscosity model, the only 
unknown is the subgrid viscosity itself which implies a closure of the form 

i^sgs = C(u). (22) 

3.3.1 Mixed scale model 

In the sequel, we focus on a subgrid-viscosity model proposed by Sagaut [15] 
having a triple dependency on the large and small structures of the resolved 
field, and the filter cutoff length. With respect to the Smagorinsky model used 
by Winckelmans et al. [13], the model proposed by Sagaut offers the advantage 
of automatically vanishing if subgrid scales are absent of the solution. This 
model, which makes up the one-parameter mixed scale family, is derived by 
taking a weighted geometric average of the models based on large scales and 
those based on the energy at cutoff. The closure is given by 

z/sgs = C,|^(u)r(gJ^A'+\ (23) 

where and 7 are the subgrid-viscosity and mixed-scale constants, is the 
resolved kinetic energy at cutoff and 

J^(u) = S(u) = S or J^(u) = V X u = aJ. (24) 
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The resolved kinetic energy at cutoff can be evaluated using the formula 

_ 1 

Qc = ^Uc,iUc,i, (25) 

where the cutoff velocity field Uc represents the high-frequency part of the 
resolved field, defined using a second filter, referred to as test filter, designated 
by the tilde symbol and associated with the cutoff length A > A 

Uc = u — u. (26) 

We note that for 7 G [0, 1], the subgrid viscosity is always defined. The con- 
stant can be evaluated by theories of turbulence in the case of statistically 
homogeneous and isotropic turbulent flow 

C, = Cl-^Cl\ (27) 

where the Smagorinsky constant — 0.18 and Cq ~ 0.20. 



3.3.2 Dynamic evaluation of the subgrid-viscosity constant 

Theoretical values of the subgrid-viscosity constant cannot be used in our case 
because they are derived if the model is used without the ADM structural con- 
tribution, that is to model the whole subgrid tensor. In order to overcome this 
issue, we introduce a dynamic procedure of Germano-Lilly type to evaluate 
this parameter as a function of space and time. Such procedure completes 
the definition of the subgrid model based on the coupling of ADM with the 
dynamic mixed scale (DMS) model, called ADM-DMS in the sequel. This 
requires the introduction of the twice-filtered Navier-Stokes equations. Ap- 
plying the test filter T^, represented by a tilde, to the filtered Navier-Stokes 
momentum equation ([9]) gives 

^ + f (n) = [f , T^] (u) + T ^ [f , g^] (u) , (28) 

which can be recast in the form 

^ + f(n) = -V-(L + r) = -V-T, (29) 

where T = L + r is an expression of the Leibniz identity referred to as multi- 
plicative Germano identity in the LES framework [16]. The components of r 
are given in Eq. and those of L by 

Lij = UiUj — UiUj, (30) 
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leading to the following expression for the subgrid tensor T corresponding to 
the twice-filtered Navier-Stokes equations 

Tij = i^Uj — UiUj. (31) 

The tensors corresponding to filtered and twice-filtered equations are modeled 
by blending ADM with the mixed scale model previously introduced. Assum- 
ing each subgrid tensor can be modeled using the same dynamic parameter 
Cd replacing the constant C^, which relies on the scale similarity hypothesis 
between test filter and primary filter cutoff lengths A and A, we obtain 



= <n* - u* u* + C^(3,,, = -2A^" '|^(u)r(gj — (32) 



= ulu* - uX u* + Cdftij, a,, = -2A |J^(u)P(gJ — 5^,-, (33) 
where cx and /3 are the subgrid- viscosity terms deprived of their constant. The 
parameter is evaluated in order to minimize the residual 

Ei^ = Lij - L™, (34) 

where L"^ = T'" - t"^. Using Eqs. (132])-([33D, Eq. §^ reads 



Eij = Lij - [{u*u* - u* u* + Cdaij) - {u*u* - u* u* + CdPij)]. (35) 

Assuming is constant over an interval at least equal to the test-filter cutoff 
length such that CdPij = CdPij, we have 

Eij = Lij - {Hij + Cdrriij), (36) 

where 

rriij = ttij — Pij and Hij = u* u* — u* u*, (37) 

which consists in a system of six independent equations leading to six possible 
different values of the constant. In a similar framework and in order to ob- 
tain a single value, Lilly [17] proposed an evaluation based on a least-squares 
minimization of the form 

^ = 0. (38) 
leading to the solution of the following single scalar equation 

Cd = (39) 

TTlij Tflij 

Smaller values than theoretical ones are expected for Cd using the previous 
dynamic procedure because of the small difference between the tensors L and 
H, only induced by the deconvolution error. Indeed, the tensor H can be 
explicitly written as 



Hij = (QN-^Ui) (Qat^Mj) - (Qjv^Mi) (Qat^Mj), (40) 
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and if the deconvolution order N oc, corresponding to Qn ^ G ^ if the 
series (fT^ is convergent, one has 

lim H^j = Lij, (41) 

which imphes that the subgrid-viscosity term vanishes if exact deconvolution 
is performed up to the grid level. This behavior of the eddy-viscosity part of 
our model, when the deconvolution order tends to infinity is strictly equivalent 
to the one observed by Sagaut et al. [31] and Stolz et al. [2] using high-pass 
filtered subgrid-viscosity models. Furthermore the relaxation term introduced 
by Stolz et al. [3,32] to stabilize their ADM-based LES has the same behavior 
in the infinite deconvolution order limit. The choice of the deconvolution order 
can be interpreted as a way of tuning the relative part taken by the subgrid- 
viscosity term which compensates the deconvolution error to minimize the 
difference between L and L™, in a least-squares sense. In the limit of going 
to infinity, the modeled subgrid stress tensor defined in Eq. (!32|) reduces solely 
to its ADM contribution 

T^l = UiUj - UiUj, (42) 
where u = P t^t u is the grid-filtered velocity. 

3.4 Particular cases of ADM-DMS 

In this section we highlight two particular cases of ADM-DMS. The first one 
is the model proposed by Zang et al. [22], based on Bardina and Smagorin- 
sky models with a dynamic evaluation of the subgrid-viscosity constant. The 
ADM-DMS expression of the subgrid tensor given by (132!) can be explicitly 
written as 

= {QN *Ui){QN*Uj) - (Qn^u,) (Qn^uj) + CdPij . (43) 

Then choosing 7 = 1, A^ = and J^(u) = S leads to 

= - % % + CdPij, Pij = -2A^\S\Sij, (44) 

which is the expression of the one-parameter dynamic mixed model. For A^ = 
0, the tensor H explicitly expressed by Eq. (l40l) reads 

Hij = UiUj — UiUj. (45) 

The second particular case of ADM-DMS is DMS, a dynamic version of the 
mixed scale model proposed by Sagaut [15]. This model is formally obtained 
by imposing Qat = in the developments of Sect. 13.3.21 which leads to H = 
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and to the following expression of the subgrid tensor 

r^^ = CdA„ A, = -2A'+^|^(u)r(gJ^5,„ (46) 
with the dynamic parameter of DMS given by 

Cd = hl^. (47) 

Without the ADM contribution, higher values of the dynamic parameter are 
expected since the difference between L and H occurring in Eq. ( l39l) disappears 
in Eq. (H7|) . This phenomenon is in direct relation with the fact that the subgrid 
viscosity term is used to model the whole subgrid tensor in this particular case. 



4 Filtering 



Filtering techniques suited to SEM and LES must preserve C°-continuity of 
the filtered variables across spectral elements and be applicable at the ele- 
ment level. In the sequel, we present a filter satisfying these constraints which 
is based on spectral techniques ensuring the element-level filtering property. 
The filtering operation is performed by applying a given transfer function in 
a modal basis. Depending on this transfer function, this filter may not be 
projective, therefore ensuring its invertibility which is a key feature needed by 
the deconvolution procedure. Hence, we will focus on the choice of the transfer 
function to fulfill this constraint. 



Description of the filter 



The modal basis introduced in the p-version of finite elements and first used 
by Boyd [18] as filtering technique is presented in its one- dimensional version, 
the extension to three dimensions being straightforward by tensor product. It 
is built up on the reference parent element VL = [—1, 1] of the SEM as 

^ ^ (48) 

0, = %(0- V2(0, 2<j<p, 

where Lj is the Legendre polynomial of degree j. Unlike the Lagrange-Legendre 
nodal basis used in our spectral element calculations, this modal basis forms a 
hierarchical set of polynomials allowing to define in an explicit and straightfor- 
ward manner a low-pass filtering procedure. Any variable v can be expressed 
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in this basis by the relation 

^(0 = E^A(0, ^efi, (49) 

j=0 

which in matrix notation reads 

V = $v, (50) 

where 

•f^.— ^ite). (51) 

The filtering operation is performed in the spectral modal space through a 
diagonal matrix K whose components are chosen in order to fulfill the required 
properties of the filter. The filtering process for a one-dimensional problem is 
expressed by 

v = $K$^v = Gv. (52) 



4 ■ 2 Transfer function 



C^-continuity, conservation of constants, invertibility and low-pass filtering 
are obtained by properly choosing the transfer function represented by the 
diagonal transfer matrix K. Imposing all these requirements to the filter could 
seem like an intractable issue but appears feasible when visualizing the modal 
basis functions presented in [33] and reported in Fig. [1] As the filter acts in 
another basis than the one used for our spectral element calculations, C°- 
continuity is preserved if the boundaries of the elements are not affected by 
the filtering procedure. One can notice that the only shape functions having 
non-zero values at the element boundaries are 0o and 0i, while (pj, j > 2 are 
bubble functions. The functions 0o and 0i are responsible for imposing the 
non-zero values on element edges. Therefore, the transfer function coefficients 
must satisfy 

Kij = Sij, i, j <l, (53) 
with Sij the Kronecker operator. If Kij verifies (l53l) . the constants are conserved 
after filtering because they are expressed as a linear combination of 0o and 0i. 
The modal filter is not projective if all diagonal coefficients Ka are non-zeros. 
The last required property is to perform low-pass filtering in frequency. As this 
modal basis forms a hierarchical set of polynomials, low-pass filtering is done 
by damping the high-degree polynomial contributions. The transfer matrix is 
expressed by 

Kij = 5ij}C{i), (54) 
with the continuous transfer function 



1+0 
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where t] is referred to as filtering rate (Fig. [2]). The transfer function is such 
that the filter verifies all the required properties previously described. The 
cutoff frequency k is arbitrarily defined by IC{k) = 1/2. Such filtering technique 
has already been used by Blackburn and Schmidt for the LES of channel flow 
using SEM [33]. In the present work, the transfer function given by Eq. fl55l) 
and depicted on Fig. [2] ensures the invertibility of the filter contrary to [33]. 
Moreover, the shape of the transfer function in Fig. [2] is similar to the one 
classically used by Stolz et al. [3,32]. However, in the SEM framework the 
transfer function is defined element by element in the spectral modal space 
which prevents from a direct comparison with the discrete filter implemented 
by Stolz et al. in [3,32]. 




Fig. 1. Bases associated with the filtering operation, shown for polynomial order 
p = 5 on the reference parent element (7 = [— 1, 1]. The Legendre polynomials (left 
column), the modal polynomials (right column). 



4- .3 Filter cutoff length 

The subgrid-viscosity term of ADM-DMS makes explicitly use of the filter 
cutoff length which needs to be defined. For a one-dimensional problem, e.g. 
in the x-direction, using the SEM, a common choice of filter width [20,21] is 

A, = (56) 

Px,c 
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Fig. 2. Transfer functions associated with Q-k for different values of the filtering rate 
r] = 1.25, 2.50, 5.00, 10.0 (left) and, H* for different values of the deconvolution 
order = 0, 1, 3, 5 with r] = 1.25 (right). Polynomial order p = 8. 



where A^. is the element size and px,c the highest polynomial degree in the 
spectral decomposition (H9|) that is the closest to the cutoff frequency k 

Px,c = max(i), such that i < k, i = 0, . . . ,p. (57) 

We notice that the filter length decreases when the element is refined and the 
polynomial degree augmented. The straightforward three-dimensional exten- 
sion for problems with rectilinear spectral elements is 



^ ^ ^ 1/3 

A{x, y, z) = {Mx)Ay{y)A,{z)Y/' ^ fA^A^AA 

\Px,c Py,c Pz,c / 



4-4 Filtering operators related to ADM 



The filtering operators Qn-^ and Hat* are defined with respect to Q-k, see Eq. 
f|T^ and f|T6|) respectively and exphcitly depend on the deconvolution order 
A^. By representing in Figure [2]the transfer function associated with T-C^-k, one 
can observe that the deconvolution error is important at the end of the modal 
spectrum, so that Hn* constitutes a high-pass filter. When increasing the 
deconvolution order N, the transfer function associated with the filter Ti^-k 
diminishes, showing the increasing accuracy of the approximate deconvolution 
procedure. 
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5 LES of the lid-driven cubical cavity flow 

5.1 General considerations 

The different LES presented hereafter refer to the flow in a lid-driven cubical 
cavity performed at Reynolds number of 12'000. The flow domain f2 consists 
in a cubical cavity such that O = (— /i, /i)^, the axis origin being assigned at 
the center of the cavity (Fig. [3]). The flow is driven by imposing a non-zero 
velocity parallel to the x-axis on the "top" wall. On the other walls, no-slip 
conditions are imposed. The moving wall will be referred to as the lid while 
the faces normal to the z-axis will be referred to as side walls. The upstream 
and downstream walls are normal to the x-axis and characterized by their 
relative position with respect to the lid motion. The remaining face parallel 
to the lid is called bottom wall. As far as the velocity imposed on the lid is 
concerned, the unit velocity induces severe discontinuities along the top edges. 
In order to avoid these defects, the imposed velocity on the lid is given by the 
polynomial expression 

u,{x, h, z) = f/o(l - {xlhYf{l - {z/hrf, Uy = u, = 0, (59) 

where Uq is a constant. The Reynolds number is defined using the maximum 
velocity Uq 




Fig. 3. Lid-driven cubical cavity. Geometry and definitions. 

Although the geometry is very simple, the flow presents complex physical 
phenomena [19,21], no direction of homogeneity and a large variety of flow 
conditions. For such Reynolds numbers, the flow over most of the domain 
is laminar and turbulence develops near the cavity walls. Its main feature is 
a large scale recirculation which spans the cavity in the 2;-direction. Aside 
this large flow structure, the relatively high momentum fluid near the lid is 
deviated by the downstream wall into a down flowing nonparallel wall jet which 
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separates ahead of the bottom wall. A region of high pressure and dissipation 
located at the top of downstream wall results from this deviation. The energy 
resulting from the impingement of the separated layer against the bottom wall 
is lost to turbulence and partly recovered by an emerging wall jet near the 
upstream wall where the flow slows down and relaminarizes during the fluid 
rise. The flow is also characterized by multiple counter-rotating recirculating 
regions at the corners and edges of the cavity. 





Time-step 


Lid vel. 


Int. time 


Nb. elements 


Polynomial degree 






n 




{Ex, Ey, Ez) 


{Px,Py,Pz) 


DNS 


0.0025 


18 


rooo 


(1,1,1) 


(128,128,128) 


LES 


0.0020 


18 


80 


(8,8,8) 


(8,8,8) 


Table 1 



Numerical and physical parameters of the DNS [19] and LES. 



The physical and numerical parameters of the DNS and the LES are gathered 
in Table [TJ The DNS constitutes the reference solution and was obtained with 
a Chebyshev collocation method on a grid composed of 129 collocation points 
in each spatial direction [19]. For LES, the spectral elements are unevenly 
distributed (Fig. Hj) in order to resolve the boundary layers along the lid and 
the downstream wall. The spatial discretization has = Ey = = 8 
elements in the three space directions with Px = Py = Pz = 8 polynomial 
degree, equivalent to 65^ grid points in total. The mesh used for LES has 
therefore twice less points per space direction than the DNS grid of Leriche 
and Gavrilakis but it is important to note that to achieve a DNS using the SEM 
would require more than 129^ grid points due to the lower order of the SEM as 
compared to the Chebyshev collocation method [23]. The space discretization 
is strictly equivalent to the one used for the LES reported by Zang et al. [22] 
for a lower Reynolds number of 10 '000. One should notice that the time-step 
for LES is slightly smaller than for the DNS which is due to different CFL 
constraints for the two different numerical schemes used, namely SEM and 
Chebyshev collocation. 

The mixed scales constant is set to 7 = 0.5 in order to have the triple de- 
pendency on the large and small structures of the resolved field as a function 
of the filter cutoff length. Furthermore, the ratio between both filtering rates 
rjr and rjg in Eq. fl55|) . corresponding to the test and primary filters Tic and 
Q-k respectively, is taken equal to two leading to a ratio of the filter cutoff 
lengths A/ A = 7/4. The parameters chosen for all LES analyzed hereafter 
are summarized in Table [21 The choice of the deconvolution order is based on 
the observations of Stolz et al. [3,32] and, Gullbrand and Chow [5] who found 
that the value N = 5 for the deconvolution order is a good compromise be- 
tween the precision in the approximate deconvolution and the computational 
cost induced by higher N in the van Cittert expansion series. This choice is 
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further justified by tlie analysis of the approximate deconvolution error devel- 
oped in Sec. 15.21 



0.5- 



> 0- 



-0.5- 



X 

Fig. 4. Spectral element grid in any plane normal to the z-direction. 

A LES based on ADM-DMS with the same parameters as the ones in Table 
121 except for JF = lJ has been carried out and has provided results extremely 
close to those from ADM-DMS with JF = S. Consequently, ADM-DMS results 
for the case = S are the only ones reported in this article. A LES based on 
DMS, see Sec. 13.41 with the same parameters as ADM-DMS for its dynamic 
mixed scale part, is also presented and compared to ADM-DMS in order to 
identify the improvement induced by coupling ADM with DMS. 

LES model N r]g rjj- ^ T 
ADM-DMS 5 1.25 2.50 0.50 S 
DMS - 1.25 2.50 0.50 S 

Table 2 

Models parameters for both ADM-DMS and DMS. 

The different LES are all started from the same initial condition, namely an 
instantaneous velocity field obtained from the DNS in the statistically-steady 
range and re-interpolated onto the spectral element grid. The projective filter 
due to this re- interpolation induces the unrecoverable loss of the subgrid scales. 

In order to verify that our mesh is coarse enough and does not resolve all 
scales of the fiow, a DNS of the lid-driven cubical cavity fiow was performed 
with SEM and with exactly the same physical and numerical parameters as 
the ones reported in Table [H One can observe on Fig. [5] that this under- 
resolved DNS (UDNS) is totally inoperative in the particular context of this 
simulation. Even first-order statistics are far from being well predicted, not 
to mention second-order ones. These results allow us to confirm the sufficient 
under-resolution of the fiow using the 65^ SEM mesh. 

We assume that a statistically-steady state is reached and time averaging will 
be taken as ensemble averaging. For any variable f , the Reynolds statistical 
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Fig. 5. In the mid-plane z/h = 0, UDNS (left column) and DNS (right column). 
Top row: contours of the x-component of the average resolved velocity field from 
—0.2 to 1 by increments of 0.01. Bottom row: contours of the y-component of the 
average resolved velocity field from —0.7 to 0.1 by increments of 0.01. Color scale 
from blue to red. Dashed contours correspond to negative levels. Levels in Uq units. 

decomposition 

v={v)+v° (61) 

introduces the time-averaged value denoted into brackets (f) and its fluc- 
tuating part v°. It is noteworthy reminding here the difference between the 
filter splitting v = v + v' and the Reynolds decomposition. As the initial con- 
dition of all LES is the same DNS instantaneous velocity field taken from 
the statistically-steady-state range, it is reasonable to also assume that LES 
will reach a statistically-steady state very quickly, if subgrid modeling is effi- 
cient [21]. These assumptions are easily verified by evaluating the total kinetic 
energy of the resolved field 

Q(u) = \ j UiUi df], (62) 

which is expected to evolve within a relatively small fiuctuation range. For 
figure O the results reported for ADM-DMS correspond to a longer dynamic 
range of 200 h/UQ time units. However, all the statistical results presented 
hereafter for both ADM-DMS and DMS are hmited to the first 80 h/Uo time 
units. The time histories of Q(u) presented in Fig. [6] for ADM-DMS and 
DMS show an evolution within the same fiuctuation range as the DNS and 
around the average value of the total kinetic energy predicted by the DNS. As 
reported by Bouffanais et al. [21] using a dynamic Smagorinsky model, which 
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is a particular case of the present DMS over 800h/Uo time units, further 
confirms the evolution of Q(u) for DMS in the long run. 



Additionally the time histories of the kinetic energy of the fluctuating resolved 
field Q{u°) presented on Fig. [7] is also evolving in the same fluctuation range 
as the DNS. The results on Fig. Eland Fig.[7]for both Q{n) and Q(u°) allow to 
conclude that the turbulent flow reaches a statistically-steady state extremely 
quickly. No transient can be clearly identified in this case. The same conclusion 
is made by Bouffanais et al. [21] for LES using more classical subgrid models. 




100 150 200 

t in h/Uo units 

Fig. 6. Total resolved kinetic energy Q(\i) in C/g/i'^ units with respect to the time t 
in h/Uo units and, for the DNS, ADM-DMS and DMS (limited to 80 time units). 
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Fig. 7. Fluctuating resolved kinetic energy Q(u°) in [/q/i^ units with respect to the 
time t in hjVQ units and, for the DNS, ADM-DMS and DMS (limited to 80 time 
units) . 

Since the scale separation used for LES leads to the removal of subgrid scales 
mainly responsible for the energy dissipation, the subgrid model has to take 
into account this phenomenon. The flow in the cavity is confined and recircu- 
lating so that the same fluid is conserved inside the cavity. Moreover, kinetic 
energy is constantly provided to it by viscous diffusion. Hence, integral energy 



260 



quantities over the flow domain such as Q{u) and Q(u°) are a direct indica- 
tion of any under- or over-dissipative character of the subgrid model, keeping 
in mind the very low numerical dissipation and dispersion of the SEM. The 
results obtained for Q{u) using ADM-DMS and DMS clearly show that the 
energy balance is achieved when using these models in this context. 



5.2 Validation of the approximate deconvolution procedure 



The flrst step towards a complete validation of the ADM-DMS model, resides 
in investigating the accuracy of the deconvolution procedure based on the van 
Cittert method, with respect to the deconvolution order N. For this purpose, 
we deflne the relative error in L^-norm between a non-flltered DNS velocity 
fleld, extracted from the DNS database of Leriche and Gavrilakis [19], and its 
deconvoluted counterpart Qn 

||U- Qiv*u||L2(j^) 

eu = TT-T: —■ (63) 

ll"llL2(fi) 



Figure [8] displays the parametric analysis of the relative error with respect to 
the deconvolution order, with the flltering rate rjg as parameter. One can no- 
tice that the van Cittert expansion series is convergent and the error increases 
with the flltering rate rjg. In practical tests the deconvolution order must be set 
lower to 30 in order to avoid having binomial coefficients of very high values 
which would inevitably induce precision errors. This numerical issue justifies 
the "apparent divergence" of the approximate deconvolution procedure for all 
filtering rate observed in Fig. [8] for large values of A^. The filtering rate is 
increased from 1 to 9 by unit increments showing that the deconvolution er- 
ror is larger with higher values of the filtering rate, which corresponds to the 
expected result. One can also notice that the error growth in the "apparent di- 
vergence" occurs earlier with lower filtering rates. It is very interesting to note 
that the error analysis is being performed using a velocity field corresponding 
to a turbulent fiow including laminar regions. The resulting deconvolution er- 
ror is clearly higher than the one obtained with a smooth analytically-defined 
field. 



5.3 A posteriori validations 



In this section, results of the LES are compared with the available reference 
results by analyzing first- and second-order statistics. The measurements re- 
ported by Leriche and Gavrilakis [19] were taken in the mid-plane z/h = 0, 
which is the statistical symmetry plane of the flow domain. For the sets of DNS 
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Fig. 8. Parametric analysis of the deconvolution error Cu with respect to the decon- 
volution order A^. The filtering rate t] is increased from 1 to 9 with unit increments. 

data, the total velocity field is considered whereas in the case of LES, only its 
resolved part is presented. In consequence, the statistical moments computed 
from the resolved field cannot be equal to those computed from the DNS. One 
solution to overcome this issue would have been to apply the same filtering 
as is used for the scale separation to the reference solution [34]. We refer the 
reader to the monograph by Sagaut [28] for more details. The statistics for 
all LES and UDNS are based on a sampling approximately 10 times smaller 
than the one of the original DNS, but about twice longer than the one of the 
LES of Zang et al. [22]; more precisely 400 samples are collected over 80h/Uo 
time units. The original reference DNS statistics were performed using 4'000 
samples extracted over an integration range of I'OOO h/Uo. Therefore, the LES 
statistics are not expected to be identical to the reference ones, especially the 
second-order ones. 

The comparisons with the DNS results are performed by plotting identical 
series of contour levels of the average velocity. Figures [9] displays the average 
values of the velocity field for DMS, ADM-DMS, and the DNS in the mid- 
plane of the cavity. This figure is complemented by the corresponding one- 
dimensional plots presented in Fig. [10] on the horizontal/vertical centerlines 
in the mid-plane z/h = 0. A rapid overview of these figures indicates that 
ADM-DMS provides results very close to the DNS references, which has to 
be compared with the UDNS results of Figure O In addition, it appears that 
ADM-DMS results are more satisfactory than those from DMS. Secondary 
corner eddies located above the bottom wall and below the lid next to the 
upstream wall are correctly resolved in the mean flow. The flow below the lid 
and near the corner with the downstream wall presents wiggles in the LES 
contours for (uy). More limited effects are noticeable for the equivalent x- 



262 



component field. We assume that these very hmited defects find their origin 
in a local too important under-resolution due to the very high shear rate near 
the downstream corner right below the lid [21]. The previous comparisons of 
ADM-DMS with the DNS and DMS for first-order moments require to be 
complemented by plotting identical series of contours of three components 
of the resolved Reynolds stress tensor. Figure [11] showcases the improvement 
achieved in terms of subgrid modeling by coupling ADM with DMS. Moreover, 
Fig. [T2HTI] provide the associated one-dimensional plots of these quantities in 
the vertical and horizontal centerlines of the mid-plane of the cavity. Indeed, 
the variations of ^l^Y^"^ and for ADM-DMS reproduce quite 

accurately the intense-fluctuations zones in the mid-plane z/h = 0, and more 
specifically in the vicinity of the downstream corner eddy. DMS appears clearly 
not as effective as ADM-DMS. The lower intensity of the Reynolds stress 
components for ADM-DMS as compared to the DNS are induced by the lower 
sampling of all LES. A longer dynamic range would produce more intense 
results as reported in [21]. 
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Fig. 9. In the mid-plane z/h = 0, DMS (left column), ADM-DMS (central column) 
and DNS (right column). Top row: contours of (ux) from —0.2 to 1 by increments 
of 0.01. Bottom row: contours of (uy) from —0.7 to 0.1 by increments of 0.01. Color 
scale from blue to red. Dashed contours correspond to negative levels. Levels in Uq 
units. 
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Fig. 10. In the mid-plane z/h = 0, DMS (dashed hnes), ADM-DMS (dotted lines) 
and DNS (solid lines). Top: {ux) on the horizontal centerline x/h = 0. Bottom: (uy) 
on the vertical centerline y/h = 0. 

5.4 Reynolds stresses production 



As mentioned by Leriche and Gavrilakis in [19], the largest Reynolds stresses 
production rates in the cavity are to be found in the primary elUptical jets 
parallel to the downstream wall, near the impact points just above the bot- 
tom wall. The budget equations of the resolved second-order moments {u°u°) 
governing the resolved Reynolds stresses, see [35,36], comprise a term named 
here Pij, defined by 



dxk 



dxk 



(64) 



and corresponding to the interaction of the mean flow and the Reynolds stress 
tensor. This quantity can be interpreted as responsible for the production 
of resolved Reynolds stresses and couples first- and second-order statistical 
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Fig. 11. In the mid-plane z/h = 0, DMS (left column), ADM-DMS (central column) 
and DNS (right column). Top row: contours of (m°^)^/^ from to 0.07 by increments 
of 0.005. Central row: contours of {ul"^)^^^ from to 0.15 by increments of 0.005. 
Bottom row: contours of from —0.0007 to 0.0065 by increments of 0.0002. 

Color scale from blue to red. Dashed contours correspond to negative levels. Levels 
in Uq units for and (Uy^)^/^ and in units for {u°ul). 

moments. 

In the specific case of the separated downstream-wall jet, the term P22 is the 
largest out of the set of Pij terms. After probing in the cavity, the maxima 
of the field P22 was found in the plane y/h = —0.9384 just at a very short 
distance above the bottom wall. The contours of the resolved Reynolds stress 
production term P22 in this plane are shown in Figure fT5l First, it can be noted 
that these contours are qualitatively very close to the ones obtained by Leriche 
and Gavrilakis in [19] and presents secondary structures of negative Reynolds 
stress production. The distribution of contours allow to clearly visualize the 
trace of the separated elhptical jets just before their impingement on the 
bottom wall. This separation is clearer for ADM-DMS than for DMS which 
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Fig. 12. In the mid-plane = 0, DMS (dashed hnes), ADM-DMS (dotted Hnes) 
and DNS (sohd hnes). on the vertical centerline y/h = (Top) and on the 

horizontal centerline x/h = (Bottom). 

shows once again the superiority of the coupled model. 
5. 5 Dynamic parameter for ADM-DMS 

In practical LES presented in the sequel, local negative values of the dynamic 
parameter are encountered. It was not found necessary to clip them as 
commonly done — e.g. in [22,33] — to conveniently get rid of locally destabilizing 
negative values. 

It is worth analyzing the variations of the dynamic parameter for ADM- 
DMS in the plane y/h = —0.9384 where the maximum of the resolved Reynolds 
stress production is found. As discussed in Sec. I3.3.2[ we expect that, by em- 
ploying ADM as the base model for the scale-similarity part of the subgrid 
stress tensor, the magnitude of the dynamic parameter will be reduced 
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Fig. 13. In the mid-plane = 0, DMS (dashed hnes), ADM-DMS (dotted Hnes) 
and DNS (sohd hnes). on the vertical centerline y/h = (Top) and on the 

horizontal centerline x/h = (Bottom). 



compared with that from the djTiamic mixed model and even more reduced 
compared with that from the dynamic Smagorinsky model [22]. This is con- 
firmed by our LES where three orders of magnitude separate the dynamic 
parameters for DMS and ADM-DMS. The distribution of contours of the av- 
erage dynamic parameter in Figure dn] appears clearly correlated with the 
contours of the resolved Reynolds stress production P22 in the same plane 
and presented in Figure [151 Indeed, the trace of the separated elliptical jets is 
discernibly apparent in Figure [T6l 

In addition, the maximum of P22 localized at the point Bq of coordinates 
x/h = 0.7874, y/h = -0.9384, z/h = -0.3371 (see Fig. EBl corresponds to a 
region of maximal values for the dynamic parameter. The time history of the 
local value of at the point Gq is reported in Figure [T71 and present a limited 
number of high-value peaks. Leriche and Gavrilakis in [19] and Bouffanais et 
al. [20] identified in this region of the cavity a pair of counter-rotating vortices 
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Fig. 14. In the mid-plane = 0, DMS (dashed hnes), ADM-DMS (dotted Hnes) 
and DNS (sohd hnes). (11° n°) on the vertical centerline y/h = (Top) and on the 
horizontal centerline x/h = (Bottom). 




Fig. 15. Contours of P22 from -0.04 to 0.11 by increments of 0.01. DMS (left), 
ADM-DMS (center), DNS (right). Plane y/h = -0.9384. Color scale from blue to 
red. Dashed contours correspond to negative levels. Levels in Uq /h units. 

responsible for the intermittent and intense production of Reynolds stresses. 
The presence of this coherent vortical structure seems to be detected by the 
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intense values of the dynamic parameter. 




Fig. 16. Contours of the average dynamic parameter from —0.001 to 0.001 for 
ADM-DMS. Plane y/h = —0.9384. Color scale from blue to red. Dashed contours 
correspond to negative levels. 
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Fig. 17. Time history of the local value of the dynamic parameter Cd for AD- 
M-DMS at the point Go whose coordinates are x/h = 0.7874, y/h = —0.9384, 
z/h = -0.3371. 



5.6 Subgrid activity 



As a next step, we are mostly interested in identifying the regions where 
turbulence occurs inside the cavity. For this purpose, we assume that if subgrid 
scales exist, the flow is locally turbulent and energy is exchanged between 
subgrid and resolved scales. In other words, the activity of the term modeling 
subgrid scales is a direct indication of the turbulence occurring in the cavity 
flow. A measure of subgrid activity is given by the subgrid energy transfer 
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defined by 



(65) 



This latter quantity is only relative in value because the dissipation induced 
by the fluid viscosity, denoted by e^, 

= 2vSijSij, Ey > 0, (66) 

is also responsible for an energetic action. It appears therefore legitimate to 
define and analyze the relative subgrid energy transfer 











+ 1 





(67) 



which is referred to as subgrid activity in the sequel. If it is close to zero, 
the energetic phenomena are mainly induced by the viscous effects showing 
that the fiuid is mainly laminar. Conversely, values close to the unit indicate 
a strong energetic action of the subgrid model refiecting that turbulence is 
mainly responsible of energy transfers. 

As one can see in Fig. [181 turbulence essentially occurs in the vicinity of 
the cavity walls but, as expected, very close to the walls the energetic action 
is essentially due to viscous effects. High values of subgrid activity are also 
identified at the bottom of the cavity and near upstream and downstream 
walls. One can also notice that subgrid activity is clearly reduced at the edges 
of the elements. This is a direct consequence of the nature of the filter which 
is not active at the element-boundaries. This issue cannot be avoided in this 
framework since C^-continuity of the variables across elements is essential for 
numerical stability and physical consistency reasons. 
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Fig. 18. Map of the average relative subgrid energy transfer {e^) from (blue) to 1 
(red) for ADM-DMS. Upstream wall x/h = (left), bottom wall y/h = (center) 
and mid-plane z/h = (right). 
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5. 7 Subgrid kinetic energy 



In order to complement the previous study of the subgrid activity and give 
further details about the importance of the subgrid terms in the ADM-DMS 
simulation, energetic quantities related to the subgrid scales are analyzed. 
For this purpose, we consider the filtered subgrid kinetic energy q' which is 
expressed — see [28, 36] — as the difference between the total filtered kinetic 
energy and the kinetic energy of the resolved field q = UiUi/2, 



where r™ is the modeled subgrid tensor defined in Eq. (!32l) . In order to provide 
the reader with deeper insight into the relative importance of the subgrid 
terms, we introduce the relative subgrid kinetic energy k, as the ratio between 
the subgrid kinetic energy and the kinetic energy of the resolved field 



As one can notice on Fig. [191 fhe average values of k reported in the plane 
z/h = 0.9384, have negative values mainly located at the top-left corner of the 
cavity. This shows that the subgrid model predicts backscattering, that is the 
energy transfer from subgrid to resolved scales. This region of inverse energy 
transfer corresponds to the region of intense subgrid activity as reported in 
Fig. [18] in the mid-plane z/h = 0. 

The importance of the subgrid terms already observed for the local energy 
fluxes and analyzed through the subgrid activity, is further confirmed by the 
presence of regions of intense (k). More precisely, four regions with high values 
of (k) and corresponding to the zones where the wall-jets detach from their 
respective walls [19] are easily identified in Fig. [13 

In order to highlight the zones where backscattering occurs, the iso-surface 
{q') = is plotted in Fig. [20] thereby defining the boundary between the 
backscattering and the forward-energy transfer regions. This figure shows that 
this phenomenon mainly occurs below the lid and in the down flowing jet next 
to the downstream wall. It also demonstrates the need for a complex subgrid 
model with such anisotropic flows containing various flow conditions and no 
direction of homogeneity. 



q' = -{uiUi - UiUi) = -Tii ~ - 



(68) 




(69) 
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Fig. 19. Map of the ratio between the resolved and the subgrid kinetic energy in the 
plane z/h = 0.9384. 




Fig. 20. Iso-surface of vanishing subgrid kinetic energy in the cavity. 
6 Conclusions 

LES of Newtonian incompressible fluid flows with ADM based on the van 
Cittert method using Legendre-SEM have been performed. A coupling with 
a dynamic mixed scale model was introduced. The coupling of the lid-driven 
cubical cavity flow problem at Reynolds number of 12 '000 with the SEM 
having very low numerical dissipation and dispersion appears to be a well 
suited framework to analyze the accuracy of the proposed subgrid model. 

The flltering operation is performed in a spectral modal space, generated by a 
hierarchical basis using the Legendre polynomials, through the application of 
a speciflcally designed transfer function. This transfer function is constructed 
in order to ensure continuity across elements, conservation of the constants, 
invertibility of the filter and to perform low-pass filtering. From the computa- 



272 



tion viewpoint, the filtering technique presented in this article, is the essential 
link between the SEM and ADM-based subgrid models. 

The validation of the deconvolution procedure performed using a DNS velocity 
sample, shows that the van Cittcrt method is convergent. Accounting for the 
reduced sampling and integration time, the LES performed with ADM-DMS 
show good agreement with the reference results. More precisely, first- and 
second-order statistics are in good agreement when compared to their DNS 
counterparts. Results for the Reynolds stresses production, coupling first- and 
second-order statistical moments, are also well predicted using this new model 
even with such reduced sampling. The analysis of the results obtained with 
DMS allows us to clearly identify the improvement induced by coupling ADM 
with DMS. Subgrid activity has been analyzed showing a qualitative corre- 
lation with the localization of small-scale structures in the cavity depicted 
in [21]. The importance of the subgrid kinetic energy as compared to the 
kinetic energy of the resolved field highlights the essential need for an appro- 
priate subgrid modeling. Furthermore, regions of backscatter are identified by 
ADM-DMS. 

All the presented results emphasize the efficiency of ADM-DMS when dealing 
with laminar, transitional and turbulent flow conditions such as those occur- 
ring in the lid-driven cubical cavity flow at Re = 12'000. 



Acknowledgements 

The authors would like to thank Dr. Roland von Kaenel from CFS Engineering, 
Lausanne, for insightful discussions. 

This research is being partially funded by a Swiss National Science Foundation 
Grant (No. 200020-101707), whose support is gratefully acknowledged. 

The DNS data were obtained on super computing facilities at the Swiss Na- 
tional Supercomputing Center CSCS and the LES data on Pleiades and Pleiades2 
clusters at EPFL-ISE. 



273 



References 



[I] S. Stolz, N. A. Adams, An approximate deconvolution procedure for large-eddy 
simulation, Phys. Fluids. 11 (1999) 1699-1701. 

[2] S. Stolz, P. Schlatter, L. Kleiser, High-pass filtered eddy- viscosity models 
for largc-cddy simulations of transitional and turbulent flow, Phys. Fluids. 

17 (2005) 065103. 

[3] S. Stolz, N. A. Adams, L. Kleiser, The approximate deconvolution model 
for large-eddy simulations of compressible flows and its application to shock- 
turbulent-boundary-layer interaction, Phys. Fluids. 13 (2001) 2985-3001. 

[4] P. Schlatter, S. Stolz, L. Kleiser, LES of transitional flows using the approximate 
deconvolution model, Int. J. Heat Fluid Fl. 25 (2004) 549-558. 

[5] J. Gullbrand, F. K. Chow, The effect of numerical errors and turbulence models 
in large-eddy simulations of channel flow, with and without explicit filtering, J. 
Fluid Mech. 495 (2003) 323-341. 

[6] B. Rembold, L. Kleiser, Noise prediction of a rectangular jet using large eddy 
simulation, AIAA J. 42 (2004) 1823-1831. 

[7] B. Shotorban, F. Mashayek, Modeling subgrid-scale effects on particles by 
approximate deconvolution, Phys. Fluids. 17 (2005) 081701. 

[8] R. von Kaenel, Large-eddy simulation of compressible flow using the finite- 
volume method, Ph.D. thesis, no. 15255, Swiss Federal Institute of Technology, 
ETH Ziirich (2003). 

[9] R. von Kaenel, L. Kleiser, N. A. Adams, J. B. Vos, Large eddy simulation of 
shock-turbulence interaction, AIAA J. 42 (2004) 2516-2528. 

[10] T. Iliescu, P. F. Fischer, Large eddy simulation of turbulent channel flows by 
the rational large eddy simulation model, Phys. Fluids. 15 (2003) 3036-3047. 

[II] C. D. Pruett, B. C. Thomas, C. E. Grosch, T. B. Gatski, A temporal 
approximate deconvolution model for large-eddy simulation, Phys. Fluids. 

18 (028104). 

[12] A. Dunca, Y. Epshtcyn, On the Stolz- Adams deconvolution model for the large- 
eddy simulation of turbulent flows, SIAM J. Math. Anal. 37 (2006) 1890-1902. 

[13] G. S. Winckelmans, H. Jeanmart, Assessment of some models for LES 
without/with explicit filtering, Geurts, Friedrich and Metais Eds., Direct and 
large-eddy simulation IV, Kluwer, 2001, pp. 55-66. 

[14] J. S. Smagorinsky, General circulation experiments with the primitive 
equations. I: The basic experiment. Month. Weath. Rev. 91 (1963) 99-165. 

[15] P. Sagaut, Numerical simulations of separated flows with subgrid models, Rech. 
Aero. 1 (1996) 51-63. 



274 



[16] M. Germane, U. Piomelli, P. Moin, W. H. Cabot, A dynamic subgrid-scale eddy 
viscosity model, Phys. Fluids A 3 (1991) 1760-1765. 

[17] D. K. Lilly, A proposed modification of the Germano-subgrid-scale closure 
method, Phys. Fluids A 4 (1992) 633-635. 

[18] J. P. Boyd, Two comments on filtering (artificial viscosity) for Chebyshev 
and Legendre spectral and spectral element methods: Preserving boundary 
conditions and interpretation of the filter as a diffusion, J. Comp. Phys. 143 
(1998) 283-288. 

[19] E. Leriche, S. Gavrilakis, Direct numerical simulation of the flow in the lid- 
driven cubical cavity, Phys. Fluids. 12 (2000) 1363-1376. 

[20] R. Bouffanais, M. O. Deville, P. F. Fischer, E. Leriche, D. Weill, Large-eddy 
simulation of the lid-driven cubic cavity flow by the spectral element method, 
J. Sci. Comput. 27 (2006) 151-162. 

[21] R. Bouffanais, M. O. Deville, E. Leriche, Large-eddy simulation of the flow in 
a lid-driven cubical cavity, Phys. Fluids. Submitted. 

[22] Y. Zang, R. L. Street, J. R. Koseff, A dynamic mixed subgrid-scale model and 
its application to turbulent recirculating flows, Phys. Fluids A 5 (1993) 3186- 
3193. 

[23] M. O. Deville, P. F. Fischer, E. H. Mund, High-order methods for incompressible 
fluid flow, Cambridge University Press, Cambridge, 2002. 

[24] Y. Maday, A. T. Patera, E. M. R0nquist, The Pjv x ¥n-2 method for 
the approximation of the Stokes problem. Tech. Rep. 92009, Department of 
Mechanical Engineering, MIT, Cambridge, MA (1992). 

[25] Y. Maday, A. T. Patera, Spectral element methods for the incompressible 
Navier-Stokes equations, Statc-of-thc-Art Survey on Computational Mechanics, 
A. K. Noor k J. T. Oden, ASME, New- York, 1989, pp. 71-142. 

[26] J. B. Perot, An analysis of the fractional step method, J. Comp. Phys. 108 
(1993) 51-58. 

[27] W. Couzy, M. O. Deville, A fast Schur complement method for the spectral 
element discretization of the incompressible Navier-Stokes equations, J. Comp. 
Phys. 116 (1995) 135-142. 

[28] P. Sagaut, Large eddy simulation for incompressible flows: an introduction, 3rd 
Edition, Springer, Berlin, 2005. 

[29] C. G. Speziale, Galilean invariance of subgrid scale stress models in large-eddy 
simulations of turbulence, J. Fluid Mech. 156 (1985) 55-62. 

[30] J. Bardina, J. H. Ferziger, W. C. Reynolds, Improved subgrid scale models for 
large eddy simulation, AIAA J. (1980) 80-1357. 

[31] P. Sagaut, P. Comte, F. Ducros, Filtered subgrid scale models, Phys. Fluids. 12 
(2000) 233-236. 



275 



[32] S. Stolz, N. A. Adams, L. Kleiser, An approximate deconvolution model for 
large-eddy simulation with application to incompressible wall-bounded flows, 
Phys. Fluids. 13 (2001) 997-1015. 

[33] H. M. Blackburn, S. Schmidt, Spectral element filtering techniques for large 
eddy simulation with dynamic estimation, J. Comp. Phys. 186 (2003) 610-629. 

[34] C. H. Moeng, J. C. Wyngaard, Spectral analysis of large-eddy simulations of 
the convective boundary layer, J. Atmos. Sci. 45 (1988) 3573-3587. 

[35] J. Mathieu, J. Scott, An Introduction to Turbulent Flow, Cambridge University 
Press, Cambridge, 2000. 

[36] S. B. Pope, Turbulent Flows, Cambridge University Press, Cambridge, 2000. 



276 



